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ABSTRACT 

Context. Ongoing measurements of the cosmic radiation (nuclear, electronic, and 7-ray) are providing additional insight 
into cosmic-ray physics. A comprehensive picture of these data relies on an accurate determination of the transport 
and source parameters of propagation models. 

Aims. A Markov Chain Monte Carlo method is used to obtain these parameters in a diffusion model. By measuring the 
B/C ratio and radioactive cosmic-ray clocks, we calculate their probability density functions, placing special emphasis 
on the halo size L of the Galaxy and the local underdense bubble of size r^. We also derive the mean, best- fit model 
parameters and 68% confidence level for the various parameters, and the envelopes of other quantities. 
Methods. The analysis relies on the USINE code for propagation and on a Markov Chain Monte Carlo technique pre- 
viously developed by ourselves for the parameter determination. 

Results. The B/C analysis leads to a most probable diffusion slope (5 = 0.861q'q4 for diffusion, convection, and reacceler- 
ation, or (5 = 0.234l[j Qg5 for diffusion and reacceleration. As found in previous studies, the B/C best-fit model favours 
the first configuration, hence pointing to a high value for S. These results do not depend on L, and we provide simple 
functions to rescale the value of the transport parameters to any L. A combined fit on B/C and the isotopic ratios 
C"Be/^Be, ^^Al/^'^Al, ^"^Cl/Cl) leads to L = Str kpc and = 120tl° Pc for the best-fit model. This value for rn is 
consistent with direct measurements of the local interstallar medium. For the model with diffusion and reacceleration, 
L = 4^5^ kpc and ru = •it.'s^ pc (consistent with zero). We vary S, because its value is still disputed. For the model 
with Galactic winds, we find that between S = 0.2 and 0.9, L varies from 0(0) to C(2) if r;, is forced to be 0, but it 
otherwise varies from 0(0) to 0(1) (with rh ~ 100 pc for all 5 > 0.3). The results from the elemental ratios Be/B, 
Al/Mg, and Cl/Ar do not allow independent checks of this picture because these data are not precise enough. 
Conclusions. We showed the potential and usefulness of the Markov Chain Monte Carlo technique in the analysis of 
cosmic-ray measurements in diffusion models. The size of the diffusive halo depends crucially on the value of the dif- 
fusion slope S, and also on the presence/absence of the local underdensity damping effect on radioactive nuclei. More 
precise data from ongoing experiments are expected to clarify this issue. 

Key words. Methods: statistical - ISM: cosmic rays 
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Almost a century after the discovery of cosmic ra- 
diation, the number of precision instruments devoted 
to Galactic cosmic ray (GCR) measurements in the 

GeV-TeV energy range is unprecedented. The GeV 7- and"bano"on-bor"ne rhave"acqurr e"d''datI,lhIt "re main "^^^^ 

ray diffuse emission is bem^ measured by the Fermi published ( Cream. lAhn et all [2001 Tracer. lA^T^nH 

satellite (|The Fermi-LAT Collaboration! |2009D, while the [20O8I : Axic lPanov et al.ll2008rPAMELA). Anti-proton s are 

TeV diffuse emission is withm reach of ground ar- also being measured (Pamela. lAdriani et Zi MbI) and 

ra^of Ccrcnko v Telescopes (e.g. Hess, [Aharonian et al. | ^^^gg^g f^^^^^ satellite and balloon experiments 

[20061 Milagro, lAbdo et al.| [20081 ) . The high-energy spec- (Ams-02, BESS-Polar). Anti-deute ron detection should be 

achieved in a few years (Ams-02, IChoutko fc Giovacchinil 
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I2008t Gaps, iFuke et all |2008[) . A complementary view 
of cosmic-ray propagation is given by anisotropy mea- 
surements from ground experime nts of high energy (e.g. , 
the Tibet Air Shower Array s, Amenomori et all l2006t 

Eas- 



200 



Supe r-Kamiokande-I det ector, iGuillian et al.l 

TOP. lAglietta et alll2009() . This multi- messenger and multi- 
energy picture will soon be completed: neutrino detectors 
are still in development (e.g.. Icecube, KmSnct), but iden- 
tifying the sources of the GC Rs should be within reach a 
few years after data collection (jHalzen et al.l[2008l) . 

All these measurements are probes to understand- 
ing and uncovering the sources of cosmic rays, the 
mechanisms of propagation, and the interaction of CRs 
with the gas and the radiation field of the Galaxy 
(jStrong et al.l l2007l) . It is important to determine the 
propagation parameters, because their value can be 
compared to theoretical predictions for the transport 
in turbulent mat 



equation with energy gains and losses. Depending on the as- 
sumptions made about the spatial and energy dependence 
of the transport coefficients, semi-analytical or fully nu- 
merical procedures are necessary to solve this set of equa- 
tions. The solution also depends on the boundary condi- 
tions, hence on the geometry of the model for the Galaxy. 

Several diffusion mo d els ar e conside red in th e lit- 
erature (I Webber et all [19921: iBlocmen et al.l 119931: 



l2001[ iBerezhko et al.l 
20081: iFarahat et al 



Strong k Moskalenkoll998t)jones et al.ll2001 : Mauri n et al 



20031 IShibata et al.ll200ek .Evoh et al. 



2008 



Ptuskin et al. 200 



n etic fields (e.g., ICasse et aL. 
Minnie et alll2007t iTautz et al.l 



2002; 



2008L 



Yan fc Lazarianll2008l and references therein), or related to 



the source spectra predicted in acceleration models (e.g 



We use a popular two-zone 
diffusion model with minimal reacceleration, where the 
Galactic wind is constant and perpendicular to the 
Galactic plane. Th e ID and 2D version of t his model are 
discus sed, e.g., in I Jones et"an (j200lt) and iMaurin et al.l 
(I2OOII ) . For the sake of legibility, the solutions are given in 
Appendix \K\ 

Below, we reiterate the assumptions of the model, and 
describe the free parameters that we constrain in this study 
(Sect. [2:4)) . 



Marcowith eTall 1 2006': Tjchiyama et all 120071: IPlagal 120081: 
Reville et al.l[200k iRcxnold^|2008| and references therein) . ' ' ''^^sport equation 



The transport and source parameters are also related to 
Galactic astrophysics (e. g., nuclear a bundances and stellar 
nucleosynthesis— ,Silbcrberg fc Tsaol [l990: Webber 1997.) 
and to dark matter indirect detection (e.g., iDonato et al.l 

120041: iDelahave et al.ll2008bh . 

In the first paper of this series (jPutze et al.ll2009l . here- 
after Paper I), we implemented a Markov Chain Monte 
Carlo (MCMC) to estimate the probability density func- 
tion (PDF) of the transport and source parameters. This 
allowed us to constrain these parameters with a sound sta- 
tistical method, to assess the goodness of fit of the models, 
and as a by-product, to provide 68% and 95% confidence 
level (CL) envelopes for any quantity we are interested in 
(e.g., B/C ratio, anti- proton flux). In Paper I, the analysis 
was performed for the simple Leaky Box Model (LBM) to 
validate the approach. We extend the analysis for the more 
realistic diffusion model, by considering constraints set by 
radioactive nuclei. The model is the minimal reaccelera- 
tion one, with a co nstant Galactic wind perpendicular t o 
the disc plane (e.g.. lJones et al.ll200ll: iMaurin et al.ll200lt) . 
allowing for a central underdensity of gas ( of a few hun- 
dreds of pc) around the solar neighbourhood (jDonato et al.l 
I2002D . 

The paper is organised as follows. In Sect. |21 we recall 
the main ingredients of the diffusion model, in particular 
the so-called local bubble feature. We briefly describe the 
MCMC technique in Scct.|3|(the full description was given 
in Paper I). We then estimate the transport parameters in 
the ID and 2D geometry. In Sect. |4j this is performed at 
fixed L (halo size of the Galaxy), using the B/C ratio only. 
The analysis is extended in Sect. |S] by taking advantage 
of the radioactive nuclei to break the well-known degener- 
acy between the parameters Kq (normalisation of the dif- 
fusion coefficient) and L. We then present our conclusions 
in Sect. [6l 

2. Propagation model 

The set of j = 1 . . .n equations governing the prop- 
agation of n CR nuc lei in the Galaxy is described in 
iBerezinskii et al.l(|l990D . It is a generic diffusion/convection 



The differential density of the nucleus j is a function 
of the total energy E and the position r in the Galaxy. 
Assuming a steady state, the transport equation can be 
written in a compact form as 

d / ■ ■dN^\ 
^dEV^' '""'dEj ■ 

The operator C (we omit the superscript j) describes the 
diffusion K{r,E) and the convection V{r) in the Galaxy, 
but also the decay rate rrad(-£') — 1/{jto) if the nu- 
cleus is radioactive, and the destruction rate Pinci ('',£') = 
S/SA/ "'iSM('')*^crinei(i') for coUisions with the interstellar 
matter (ISM), in the form 



Cir, E) = -V- {KV) + V-V + P^ad + P;, 



(2) 



The coefficients b and c in Eq. ([T|) are respectively first 
and second order gains/losses in energy, with 



b(^r,E) = A. 



v.v 



Eh 



2m + Ek 
m + Ek 



E 



X Kr, 



c{r,E) ^ 13^ X K. 



pp- 



(3) 



(4) 



In Eq. ([3]), the i onisation and Coulomb e nergy losses 
arc taken f rom 'Mannheim fc SchlickeiseJ (|1994[ ) and 
iStrong fc Mosk alcnko (1998,1 . The divergence of the 
Galactic wind V gives rise to an energy loss term re- 
lated to the adiabatic expansion of cosmic rays. The last 
term is a first order contribution in energy from reacceler- 
ation. Equation @ corresponds to a diffusion in momen- 
tum space, leading to an energy gain. The associated dif- 
fusion coefficient Kpp (in momentum space) is taken from 
the mod el of minimal reacc eleratio n by the interstellar t ur- 
bulence (lOsborne fc Ptuski n 1988: ISeo fc Ptuskin|[T99l . It 
is related to the spatial diffusion coefficient K by 



(52) (4 -,5)' 



(5) 
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Fig. 1. Sketch of the model: sources and interactions (in- 
cluding energy losses and gains) are restricted to the thin 
disc oc 2hS{r). Diffusion K and convection Vc transport nu- 
clei in both the disc (half-height h) and the halo (half-height 
L). The Galaxy radial extension is R. The local bubble is 
featured to be a cavity of radius in the disc devoid of 
gas. 



where Va is the Alfvenic speed in the medium. 

The source term is a combination of i) 

primary sources q^{r,E) of CRs (e.g., super- 
novae), ii) secondary fragmentation- induced sources 
ET""'' riisMir)val-;;^iE)N''ir,E), and iii) secondary 

decay-induced sources -^'^(''; ^)/(7'^o^"')- partic- 
ular, the secondary contributions link one species to all 
heavier nuclei, coupling together the n equations. However, 
the matrix is triangular and one possible approach is to 
solve the equation starting from the heavier nucleus (which 
is always assumed to be a primary). 

2.2. Geometry of the Galaxy and simplifying assumptions 

The Galaxy is modelled to be a thin disc of half-thickness h, 
which contains the gas and the sources of CRs. This disc is 
embedded in a cylindrical diffusive halo of half-thickness L, 
where the gas density is assumed to be 0. CRs diffuse into 
both the disc and the halo independently of their position. 
A constant wind Vc perpendicular to the Galactic plane 
is also considered. This is summarised in Fig. [T] (see next 
section for the definition of r/j ) . 



We r ise the Sjz) approx i mation introdt iced in iJones 
1979|) . iPtuskin fc Soutoull (|l990t) . and IWebber et al 



19921) . Considering the radial extension R of the Galaxy 
to be either infinite or finite leads to the ID version or 2D 
version of the model, respectively. The corresponding sets 
of equations (and their solutions) obtained after these sim- 
plifications are presented in Appendix [X] These assump- 
tions allow for semi-analytical solutions of the problem, 
as the interactions (destruction, spallations, energy gain 
and losses) are restricted to the thin disc. The gain is in 
the computing time, which is a prerequisite for the use 
of the MCMC technique, where several tens of thousands 
of models are calculated. These semi-analytical models re- 
prod uce all salient features of fu ll numerical approaches 
(e.g.. IStrong fc Moskalenkol 119981 ). and they are useful for 
systematically studying the dependence on key parame- 
ters, or some system atics of the parameter determination 
(jMaurin et al.ll2010D . 

We note that most of the results of the paper are 
based on the ID geometry (solutions only depend on z), 
which is less time-consuming than the 2D one in terms of 



computing timcQ- The parameter degeneracy is also more 
easily e xtracted a nd understood in this case (jJones et al.l 
[2001: Maurin et al. 2006). Nevertheless, the results for the 
2D geometry are also reported, as it has been used in 
a ser i es of studies inspecting stable nuclei (iMaurin et a.1.1 
I2001L |2002[ ). /3-radi oactive nuclei (iDonato et al.l I2002D . 
standard a nti-nuclei ("Donato et al.| l200ll I2008L l2009t) and 
positrons (jDelahave et al. .2008a) . It has also been used 
to set constraints on dar k matter annihi lations in anti- 
nuclei (jDonato et al.l[2003 ). and positrons (jDelahave et al.l 
l2008b). The reader is refe rred to these papers, and espe- 
cially Maurin et al.l (|200l[ ) for more details and references 
about the 2D case. 



2.3. Radioactive species and the local bubble 

Our model does not take into account all the observed ir- 
regularities of the gas distribution, such as holes, chimneys, 
shell-like structures, and disc flaring. The main reason is 
that as far as stable nuclei are concerned, only the average 
grammage crossed is relevant when predicting their flux 
(which motivates LBM) . As such, the thin-disc approxima- 
tion is a good trade-off between having a realistic descrip- 
tion of the structure of the Galaxy and simplicity. 

However, the local distribution o f gas affects the flux 
calculation of radioact i ve species (iPtuskin et al.l Il997l: 
IPtuskin fc Soutoull [l998l : iDonato et al.ll2002f ). Weconsider 
a radioactive nucleus that diffuses in an unbound volume 
and decays with a rate 1/(7x0). In spherical coordinates, 
appropriate to describe this situation, the diffusion equa- 
tion reads 

- KAr-G + — = S{r) . (6) 



The solution for the propagator G (the flux is measured at 
r = for simplicity) is 



G(r') cx 



(7) 



Secondary radioactive species, such as ^"^Be, originate from 
the spallations of the CR protons (and He) with the ISM. 
We model the source term to be a thin gaseous disc, except 
in a circular region of radius at the origin. In the 6{z) 
approximation (see Fig. [ij and in cylindrical coordinates. 



Q{r, z) cx 9(r - ru)5{z) 



(8) 



where is the Heaviside function. The flux of a radioac- 
tive species is thus given by (we rewrite the propagator in 
cylindrical coordinates) 



iV(r = z = 0) cx 



+ CXD 



G( v/r'2 +z'2) Q(r', z') r'dr'dz'. (9) 



The ratio of the flux calculated for a cavity/hole rh to that 
of the flux without hole {r^ = 0) is 



Nr 



exp 



rh=0 



-rh 



exp 



7^ ) • (10) 



The 2D solution is based on a Bessel expan- 
sion/resummation [see Eq. (IA.14|) ]. For each Bessel order, 
an equation similar to that for the ID geometry needs to be 
solved. Nine Bessel orders are in many cases enough to ensure 
convergence (jMaurin et al.l [20011 ). but at least 100 orders are 
required in the general case, which multiply the computing 
time by roughly the same amount. 
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The quantity ^i-ad = \JK^t^ is the typical distance on 
which a radioactive nucleus diffuses before decaying. Using 
K sa 10^* cm^ s~^ and t w 1 Myr, the diffusion length is 
'rad ~ 200 pc. Hence, in principal, any undcrdcnsity on a 
scale Th ~ 100 pc about the Sun leads to an exponential 
attenuation of the flux of radioactive nuclei. This attenu- 
ation is both energy-dependent and species-dependent. It 
is energy-dependent because it decreases with the energy 
as the time-of-flight of a radioactive nucleus is boosted by 
both its Lorentz factor and the increase in the diffusion 
coefficient. It is species-dependent because nuclei half-lives 
for the standard Z < 30 cosmic-ray clocks range from 0.307 
Myr for ^Sci to 1.51 Myr for lORe. 

In this paper, we model the local bubble to be this sim- 
ple hole in the gaseous disc, as shown in Fig. [TJ The expo- 
nential decrease in the flux of this modified DM, as given 
by Eq. ([TU)) . is directly plugged into the solutions for the 
standard DM (rh = 0). In principle, i) the hole has also 
an impact on stable species as it decreases the amount of 
matter available for spallations, and ii) in the 2D geometry, 
a hole at Rq = 8 kpc bre aks down the cylind rical geome- 
try. However, in practice, iDonato et al.l (l2002[ ) found that 
the first effect is minor, and that the hole can always be 
taken to be the origin of the Galaxy (the impact of the R 
boundary being negligible for radioactive species). 

Other sub l eties exist, which were not considered in 
iDonato et al.l (|2002[ ). Indeed, the damping in the solar 
neighbourhood — combined with the production of the ra- 
dioactive species matching the data at low energy — means 
that at intermediate GeV/n energies, the flux of this ra- 
dioactive species is higher in the modified model (with 
fh 7^ 0) than in the standard one (with — 0). It also 
means that everywhere else in the Galactic disc, at all en- 
ergies, the radioactive fiuxes are higher in the modified 
model (with damping). There are two consequences: i) all 
spallative products from these radioacti ve nuclei originate 
in an effective diffusion region in the disc (jTaillet fc MaurinI 
[200l . the size of which may be much larger than the size of 
the underdense bubble. In this case, these products ought 
to be calculated from the undamped fluxes; ii) the decay 
products of these radioactive nuclei (e.g., ^^B, which orig- 
inates from the /3-decay of ^°Be) are stable species that 
originate in an effective diffusion sphere (decay can occur 
not only in the disc, but in the halo). Both these effects 
must be considered because their contributions potentially 
affect the calculation, e.g., of the B/C and Be/B ratios (by 
means of the B flux), which are used to fit the models. We 
confirm that taking spallative products from the damped 
or undamped radioactive fluxes left these ratios unchanged. 
On the other hand, for the decay products, the effect is of 
the order of 1 — 10%, which is in general enough to change 
the values of the best-fit parameters. However, the average 
fiux (over the effective diffusion zone) from which the decay 
products originate lies between the damped and undamped 
values: the lower the effective diffusive sphere, the closer 
the flux is to the damped one. In particular, at low energy, 
when convection is allow ed, the diffusion zone can be small 
(|Taillet fc Maurinll2003D . 

To keep the approach simple, we use here the damped 
flux of radioactive species for all s pallative and decay prod- 
ucts (as was implicitly assumed in lDonato et al.ll2002|) . This 
approach is expected to provide the maximal possible size 
for r/i (if a non-null value is preferred by the fit). 



2.4. Input ingredients and free parameters of the study 

2.4.1. Gas density 

The gas density scale height strongly varies with r de- 
pending on the form considered — neutral, molecular, or 
ionised (see, e.g., .Fc rriere 200TI) . We use the surface den- 
sity measured in the solar neighbourhood as a good es- 
timate of the average gas in the Galactic disc. We set 
nisM ~ 1 cm^'^, which corresponds to a surface d ensity 
SisM = 2/inisM - 6 X 10^" cm'^ (jFerrierd [200l . The 
number fraction of H and He is taken to be 90% and 
10%, respectively. The ionised- hydrogen space- averaged 
density may be identified with the free-electron space- 
averaged density, which is the sum of the contributions 
of H ii regions and the diffuse component (|G6mez et al.l 
120011 : lFerrierel[200l . The intens ity of the latter is well 
measured 0.018 ± 0.002 cm~ 3 (jBerkhuiisen et all 120061: 
iBerkhuiisen fc Miilleil |2008[) . whereas t he former depends 
stron gly on the Galactocentric radius r ([Anderson fc Banial 
120091 ). For the total electron density, we choose to set 
(n^-) = 0.033 and - lO^K (iNordgren et al.l[T99l . 

The disc half-height is set to he h = 100 pc. It is not 
a physical parameter per se in the S{z) approximation, al- 
though it is related to the phenomena occurring in the thin 
disc. Physical parameters are related to the surface den- 
sity, which is easily rescaled from that calculated setting 
h — 100 pc (should we use a different h value). In the 2D 
geometry, the boundary is set to be i? = 20 kpc and the 
sun is located at Rq = 8.0 kpc. 

2.4.2. Fragmentation cross-sections 

In Paper I, the sets of fragmentation cross-sections 
were taken from the semi-e mpirical formulatio n of 
IWebber et all (ll990D up dated in lWebber et al.l ([1998D (see 
also Maurin et al.ll200ll and references therein) . In this pa- 
per, they are replaced by the 2003 version, as given in 
IWebber et all ([2003D . S pallations on He are calculated with 
the parameterisation of IFerrando et all i 



2.4.3. Source spectrum 

We assume that a universal source spectrum for all nu- 
clei exists, and that it has a simple power-law description. 
As in Paper I, we assume that Q{E) oc [i^R~°'. The pa- 
rameter a is the spectral index of the sources and 77 en- 
codes the behaviour of the spectrum at low energy. The 
normalisations of the spectra arc given by the source abun- 
dances qj, which are renormalised during the propagation 
step to match the data at a specified kinetic energy per nu- 
cleon (usually ~ 10 GeV/n). The correlations between the 
source and the transport parameters and their impact on 
the transport parameter determination were discussed in 
Paper I. In this s tudy, we set 77 = —1 and ^ ^ a + 5 — 2.65 
(lAve et al.ll2008D . Constraints on the source spectra from 
the study of the measured primary fluxes are left to a sub- 
sequent paper (Donato et al., in preparation). 



2.4.4. Free parameters 

We have two geometrical free parameters 

— L, the halo size of the Galaxy (kpc); 
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— rji, the size of the local bubble (kpc), which is most of 
the time set to be (to compare with models in the 
literature that do not consider any local underdensity) ; 

and four transport ones 

— Kq, the normalisation of the diffusion coefficient (in unit 
of kpc^ Myr^i); 

— 5, the slope of the diffusion coefficient; 

— Vc, the constant convective wind perpendicular to the 
disc (km s~^); 

— Va, the Alfvenic speed (km s~^) regulating the reaccel- 
eration strength [see Eq. ([S])]. 



Table 1. Classes of models tested in the paper. 



The diffusion coefficient is taken to be 



(11) 



3. MCMC 



The MCMC method, based on the Bayesian statistics, is 
used here to estimate the full distribution — the so-called 
conditional probability-density function (PDF) — given 
some experimental data (and some prior density for these 
parameters). We summarise below the salient features of 
the MCMC technique. A detailed description of the method 
can be found in Paper I. The issue of the efficiency, which 
was not raised in Paper I, is discussed in Appendix [C] 

The Bayesian approach aims to assess the extent to 
which an experimental dataset improves our knowledge of 
a given theoretical model. Considering a model depending 
on m parameters 



(12) 



we wish to determine the PDF of the parameters given the 
data, P(0|data). This so-called posierior probability quan- 
tifies the change in the degree of belief one can have in the 
771 parameters of the model in the light of the data. Applied 
to the parameter inference, Bayes theorem reads 



P(0|data) 



P(data|6/) • P{e) 
P(data) 



where P(data) is the data probability (the latter does not 
depend on the parameters and hence, can be treated as 
a normalisation factor). This theorem links the posterior 
probability to the likelihood of the data C{6) = P(data|0) 
and the so-called prior probability, P{0), which indicates 
the degree of belief one has before observing the data. 
The technically difficult point of Bayesian parameter es- 
timates lies in the determination of the individual poste- 
rior PDF, which requires an (high-dimensional) integration 
of the overall posterior density. Thus an efficient sampling 
method for the posterior PDF is mandatory. 

In general, MCMC methods attempt to studying 
any target distribution of a vector of parameters, here 
P(0|data), by generating a sequence of n points/steps 
(hereafter a chain) 



{e^}^ 



{Oi, 02, . . . , On}- 



(14) 



Each 9i is a vector of m components, e.g., as defined in 
Eq. ([T2|). In addition, the chain is Markovian in the sense 
that the distribution of 9n+i is influenced entirely by the 
value of On- MCMC algorithms are developed to ensure 



Model 


Transport parameters 


Description 


I 


{Ao, S, Vc} 


Diffusion + convection 


II 


{Ko, 5, Va} 


Diffusion + reacceleration 


III 


{Ko, S, Va} 


Diff. + conv. + reac. 



Table 2. Most probable values for B/C data only (L 
4 kpc). 



Model 7^0 X 10^ 
Data (kpc^Myr 



Vc Va 

(kms~^) (kms~^) 



I-F 


n 42+0-"^ 

^■^^-0.04 


qQ+O-02 

U.30_Q Q3 


i3.5«:^3 




II-F 


n 7+0-3 

^- ' -0.2 






73tl 


III-F 








38tl 



that the time spent by the Markov chain in a region of the 
parameter space is proportional to the target PDF value 
in this region. Here, the prescription used to generate the 
Markov chains is the so-called Metropolis-Hastings algo- 
rithm, which ensures that the stationary distribution of the 
chain asymptotically tends to the target PDF. 

The chain analysis is based on the selection of a subset 
of points from the chains (to obtain a reliable estimate of 
the PDF). Some steps at the beginning of the chain are 
discarded (burn- in length). By construction, each step of 
the chain is correlated with the previous steps: sets of inde- 
pendent samples are obtained by thinning the chain (over 
the correlation length). The fraction of independent sam- 
ples measuring the efficiency of the MCMC is defined to be 
the fraction of steps remaining after discarding the burn- 
in steps and thinning the chain. The final results of the 
MCMC analysis are the target PDF and all marginalised 
PDFs. They are obtained by merely counting the number 
of samples within the related region of parameter space. 



(13) 4. Results for stable species (fixed halo size L) 



For stable species, the degeneracy between the normalisa- 
tion of the diffusion coefficient Kq and the halo size of the 
Galaxy L prevents us from being able to constrain both pa- 
rameters at the same time. We choose to set L = 4 kpc (we 
also set rii ~ 0, i.e., standard DM). The free transport pa- 
rameters are {Kq, S, Vc, Va}- The classes of models consid- 
ered are summarised in Table [TJ The reference B/C dataset 
(denoted dataset F) used for the analysis is described in 
Appendix ID. II 

4.1. PDF for the transport parameters 

We begin with the PDFs of the parameters based on the 
B/C constraint (dataset F) for the various classes of models 
(I, II, or HI). The PDFs are shown in Fig. [2] 

The first important feature is that the marginal distri- 
butions of the transport parameters (diagonals) are mostly 
Gaussian. From the off-diagonal distributions, we remark 
that Kq and i5 are negatively correlated. This originates in 
the low-energy relation K{E) oc KqR^ , which should re- 
main approximately constant to reproduce the bulk of the 
data at GeV/n energy. The diffusion slope 5 is negatively 
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A. Putze et al.: An MCMC technique 




0.22 0^4 0.09 0.095 0.1 0.105 




Fig. 2. From top to bottom: posterior PDFs of models I, II, 
and III using the B/C constraint (dataset F). The diagonals 
show the ID marginalised PDFs of the indicated parame- 
ters. Off-diagonal plots show the 2D marginalised posterior 
PDFs for the parameters in the same column and same line 
respectively. The colour code corresponds to the regions of 
increasing probability (from paler to darker shade) , and the 
two contours (smoothed) delimit regions containing, respec- 
tively, 68% and 95% (inner and outer contour) of the PDF. 



Model 
Data 


K^"""- X 10^ 
(kpc^ Myr"^) 




best 

(km s'^) 


■^bcst 

(km s-^) 


xVd.o.f 


I-F 


0.42 


0.93 


13.5 




11.2 


II-F 


9.74 


0.23 




73.1 


4.68 


III-F 


0.48 


0.86 


18.8 


38.0 


1.47 



correlated with 14, which is related to a smaller 5 being 
obtained if more reacceleration is included. On the other 
hand, the positive correlation between 5 and Vc indicates 
that larger 5 are expected for larger wind velocities. 

We show in Table [H the most probable values of the 
transport parameters, as well as their uncertainties, corre- 
sponding to 68% confidence levels (CL) of the marginalised 
PDFs. The precision to which the parameters are obtained 
is excellent, ranging from a few % to 10% at most (for 
the slope of the diffusion coefficient 5 in III). This corre- 
sponds to statistical uncertainties only. These uncertainties 
are of the order of, or smaller than systematics generated 
from uncertainties in the input ingredients (see details in 
iMaurin et ai1l2010t) . 

As found in previous studies (e.g.. lLionetto et aI1l2005[ ). 
for pure diffusion/reacceleration models (II), the value of 
the diffusion slope 5 found is low (w 0.23 here). When con- 
vection is included (I and III), 5 is large (sa 0.8 — 0.9) . This 
scatter in 8 was already observed in lJones et al.l ()2001|) , who 
also studied different classes of models. The origin of this 
scatter is consistent with t he aforementioned correlations 
in the parameters (see also iMaurin et aLll2010l ). 

The best-fit model parameters (which are not always 
the most probable ones) are given in Table |31 along with 
the minimal value per degree of freedo m, y^^^/d.o.f 
(last column) . As found in previous analyses (jMaurin et al.l 
I2OOIL |2002[) . the DM with both reacceleration and convec- 
tion reproduces the B/C data more accurately than with- 
out: xVd-o.f= 1.47 for III, 4.90 for II, and 11.6 for I. The 
B/C ratio associated with these optimal values are dis- 
played with the data in Fig. [3] We note that the poor fit for 
II (compared to III) is explained by the departure of the 
model prediction from high-energy IIEAO-3 data. 

4.2. Sensitivity to the choice of the B/C dataset 

For comparison purposes, we now focus on several datasets 
for the B/C data. Low-energy data points include ACE 
data, taken du r ing t he solar minimum period 1997-1998 
(|de Nolfo et al.l l2006l ). Close to submissio n of this pa- 
per, another ACE analysis was published (iGeorge et al.l 
I2009D . The 1997-1998 data points were reanalysed and 
complemented with data taken during the solar maxi- 
mum period 2001-2003. The AMS-01 also provided B/C 
data covering almost the same range as the HEAO-3 
data (jTomassetti fc AMS-01 CoUaborationl I2009D . Hence, 
for this section only, we attempt to analyse other B/C 
datasets that include these components: 

- A: HEAO-3 [0.8 - 40 GeV/n], 14 data points; 

- C: HEAO-3+I0W energy [0.3-0.5 GeV/n], 22 data points; 

- F: HEAO-3 + low + high energy [0.2 - 2 TcV/n], 31 data 
points; 
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IMP7-S (•1>=250 MV) 
Spacelab-2 (*=200 MV) 
HEAO-3 (fc250 MV) 
Voyagerl&2 (<I>=225 MV) 
ACE 97-98 (*=225 MV) 
CREAM (■I>=425 MV) 




Table 4. Best-fit model parameters based on different B/C 
datasets. 



Model I-F (X- Mol=n.2) 
. Model II-F ix- /dof=4.68) 
. Model III-F (x^ /dof=1.47) 



10"' 1 10 10- 10- 

Ek/n (GeV/n) 

Fig. 3. Best-fit ratio for model I (blue-dotted line), II 
(red-dashed line), and model III (black-solid line) using 
dataset F: IMP7-8, Voyagerl&2, ACE-CRIS, HEAO-3, 
Spacelab, and CREAM. The curves are modulated with 
$ = 250 GV (and $ = 225 GV at low energy). The corre- 
sponding best-fit parameters are gathered in Table |3l 



.9 0.35 



U 

3 0.3 



IS 



AMS-Ol (Tomassctii e 
• HEAO-3 (Engelmann t 
< Spacelab-2 (Swordy cl 
r CREAM 04 (Aim ctal. 



n ^ 

;t al. 2009) 
et al. 1990) 
al. 1990) 
2008) 



+ 1MP7-S (Gaicia-Munoz ct ill. 1987) 
VoyLigcrl&2 (Lukasiak 1999) 

A ACE 97-98 (De noifo et al. 2006) 
ACE 97-98 (George et al. 2009) 

O ACE 01-03 (George etal. 2009) 



1 I" Ek/n(GeV/n) 

Fig. 4. B/C data used in this section. Shown are the IS 
data (rescaled from TO A data using Ej.^ = E'^'^^ -I- $, see 
Paper I). For several experiments, in addition to the error 
bars in the ratio, we display the energy interval from which 
the central energy point is obtained. 



- Gl: as F, but with new ACE 1997-1998 data, 31 data 
points; 

- G2: as F, but with new ACE 2001-2003 data only, 31 
data points; 

- Gl/2: using both 1997-1998 and 2001-2003 ACE data, 
37 data points; 

- H: as F, but HEAO-3 replaced by AMS-Ol data, 27 data 
points. 

The data are shown in Fig. 01 Thanks to the high level 
of modulation for the 2001-2003 ACE data, the IS (demod- 
ulated) B/C ratio covers nicely the gap between HEAO- 
3 and lower energy data. HEAO-3 and AMS-Ol data also 
show consistency across their whole energy range . 

The best-fit model parameters for these data are shown 
in Table |31 The low-energy data play an important part in 
the fitting procedure: S decreases by 0.1 when going from 
HI- A to III-C, and the diffusion normalisation is decreased. 
When the CREAM data at higher energy are taken into ac- 
count (III-F). the best-fit diffusion slope S again becomes 



Model -/i'o'"="* X 10^ 
Data (kpc^Myr"^) 




T r best 
''^ C 

(kms^^) 


T /best 

(kms-^) 


xVd.o.f 


III- A 2.51 


1.00 


21.7 


35.4 


2.11 


III-C 0.43 


0.89 


18.9 


36.7 


1.72 


III-F 0.48 


0.86 


18.8 


38.0 


1.47 


III-Gl 0.53 


0.84 


18.0 


37.4 


1.80 


III-G2 0.46 


0.85 


20.0 


39.6 


2.73 


III-Gl/2 0.53 


0.83 


19.0 


39.1 


2.94 


III-H 1.85 


0.51 


18.1 


54.1 


0.25 



slightly lower (from 0.89 to 0.86), but CREAM data un- 
certainty is still too important to be conclusive. The im- 
pact of the low-energy ACE reanalysed data points is seen 
when comparing III-F with III-Gl: the scatter between the 
derived best-fit parameters is already of the order of the 
statistical uncertainty (see Table [2|) . The data taken either 
during the solar minimum period (Gl) or the solar maxi- 
mum period (G2) cover a different energy range (see Fig.|3]). 
The Xmin ^'^^ ^2 is greater, which is not surprising, given 
the abnormal trend followed by these data (empty circles in 
Fig. U). Nevertheless, it is reassuring to see that they lead 
to consistent values of the transport parameters. 

If wc now replace the HEAO-3 data with the AMS- 
Ol data, the impact on the fit is striking: the best-fit dif- 
fusion slope 6 goes from 0.86 to 0.51. As discussed in 
aurin et "all ()2010() . HEAO-3 data strongly constrain the 
slope towards S « 0.8, even if there is a systematic energy 
bias in the HEAO-3 data themselves. From the AMS-Ol 
data, we see that there could be a way of reconciling the 
presence of a Galactic wind and reasonable values of S. 
However, the large error bars in AMS-Ol data, reflected by 
the low x^/d.o.f value, does not allow to draw firm conclu- 
sions. Data in the same energy range from PAMELA would 
be helpful in that respect. Moreover, high energy data from 
subsequent CREAM flights or from the TRACER experi- 
ments will be a crucial test of the diffusion slope: at TeV 
energies, diffusion alone is expected to shape the observed 
spectra, so that the ambig uity with the effect of conv ection 
or reacceleration is lifted (jCastellina &: Donatcill200"5h . 

4.3. Comparison of trends for the DM and for the LBM 

For completeness, we briefly comment on the similarities 
and differences between the results found here and in 
Paper I. To follow the organisation of the previous sec- 
tions, the comparison with the LBM is discussed for differ- 
ent classes of models (I, II, and HI), and then for different 
datasets (A, B, and C). We note that the best-fit values 
presented below differ slightly for those given in Paper I, as 
an updated set of production cross-section is used. 

We recall that in the LBM (see Paper I), the free pa- 
rameters are the normalisation of the escape length Aq, 
(5, a cut-off rigidity Rq, and a pseudo-Alfvenic speed Va- 
The latter is linked to a true speed by means of Va = 
Va X (hLf^, i.e., Va = O.4I/2 Va for h = 0.1 kpc and 
L — 4 kpc. The diffusion coefficient at 1 GV is related 
to the escape length by means of Kq ~ 0.5 c x fiL/Xo, 
where we use n = 2hnm = 1.34 x 10""^ g cm~^, leading 
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Table 5. Best-fit parameters on B/C data for the LBM. 



Model 
Data 


i^o"""' X 10^ 
(kpc^ Myr-i) 


^bcst 


Ro 
(GV) 


T r best 

(km s-i) 


x7d.o.f 


I-F 


2.36 


0.56 


5.70 




5.52 


II-F 


5.26 


0.38 




65.1 


1.78 


III-F 


4.19 


0.43 


2.94 


53.9 


1.56 


III-A 


2.32 


0.57 


4.40 


11.3 


2.71 


III-C 


4.13 


0.44 


3.10 


53.6 


2.26 


III-F 


4.19 


0.43 


2.94 


53.9 


1.56 



to Ko (kpc^ Myr-i) w 0.82/Ao (g cm-^). The LBM pa- 
rameters gathered in Table [S] are obtained from the above 
conversions, to ease the comparison with the DM results. 

For the different classes of models (I, II and III), a com- 
parison of Table |3] with the first three rows of Table [5] in- 
dicates that the same trend is found. For instance, model I 
(without reacceleration) has a larger d than those with, and 
model II (without convection/rigidity-cutoff) has a smaller 
6 than those with. The slope for model III (with both 
convection and reacceleration) is in-between. This effect is 
more marked for the DM than for the LBM. We note that 
model II (with reacceleration but without convection) is al- 
most consistent with a Kolmogorov spectrum of turbulence, 
but is inconsistent with the data. Concerning the different 
datasets (A, C, and F), again, the same trend as for the 
LBM is found (compare Table [Hand the last three rows of 
Table [5|). 

The most striking difference between the two models 
(LBM and DM) concerns their 6 values. This difference can 
be explained in terms of non-equivalent parameterization 
of th e low-energy transport coefficient (see iMaurin et al.l 
[20T0l for more details). Apart from this, both the value of 
the Alfvenic speed and the normalisation of the diffusion 
coefficient Kq in the two cases are fairly consistent when 
similar values of 6 are considered. 



4.4. Dependence of the parameters with L 

All the previous conclusions were derived for L = 4 kpc, but 
hold for any other halo size. The evolution of the transport 
parameters with L is shown in Fi g. [5] (the best-fit va lues 
are consistent with those found in lMaurin et al.l [20021) . In 
the three upper figures, we have superimposed the observed 
dependence a parametric formula. 

For Kq (top panel) , the formula can be understood if we 
consider the grammage of the DM. In the purely diffusive 
regime, we have Aosc oc L/K. This means that when we 
vary L, to keep the same grammage in the equivalent LBM, 
we need to vary Kq accordingly. We find that Kq = 1.08 x 
10-3(L/1 kpc)i °6 kpc^ Myr-i instead of Kq cx L. The 
origin of the residual L^ '^^ dependence is unclear. It may 
come from the energy loss and gain terms. 

For the reacceleration, the interpretation is also simple. 
From Eq. (O, Va should scale as y^K^, so that Va oc Vl. 
We find Va = 18.21(L/1 kpcf-^^ km s'^. This is ex- 
actly Va OC 



^Kq, with the dependence \/ L'^ og as above. 
The quantities Vc and S are roughly constant with L. The 
Xmin surface is rather flat, although a minimum is observed 
around L ^ 15 kpc (the presence of a minimum may be 
related to the presence of the decayed ^°Be into ^"^B in 
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Fig. 5. Best-fit parameters (III-F) as a function of the halo 
size of the Galaxy (blue circles). From top to bottom: Kq, 
Va, Vc, S, and the associated Xmm- the first three figures, 
a parametric function matching the observed dependence 
is shown (dashed- red line). 



the B/C ratio). This fiatness is a consequence of the de- 
generacy of Kq/L when only stable species are considered. 
Consequently, an MCMC with L as an additional free pa- 
rameter does not converge to the stationary distribution. 
A sampling of the Galactic halo size is possible if radioac- 
tive nuclei are considered to lift the above degeneracy (see 
Sect. ED. 

4.5. Summary of stable species and generalisation to the 2D 
geometry 

The transport parameters for both LBM (Paper I) and ID 
DM, when fitted to existing B/C data, are consistent with 
both convection and reacceleration. The correlations be- 
tween the various transport parameters, as calculated from 
the MCMC technique, are consistent with what is expected 
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Table 6. Best-fit model parameters on B/C data: ID versus 
2D DM (L = 4 kpc). 



Model i^o^"^' X 10" i''"^* V^"^ xVd.o.f 
Data (kpc^Myr'^) (kms~^) (km s~^) 



ID 


II-F 


9.74 


0.23 




73.1 


4.68 


2D 


II-F 


8.56 


0.24 




68.6 


4.67 


ID 


III-F 


0.48 


0.86 


18.9 


38.0 


1.47 


2D 


III-F 


0.42 


0.86 


18.7 


35.5 


1.46 



from the relationsh ips between DMs and the LBM (e.g., 
iMaurin et al.l[2006|) . From the B/C analysis point of view, 
it implies that even if we are unable to reach conclusions 
about the value of S (see aurm "eTalllMl . once this 
value is known, all other transport parameters are well con- 
strained. 

The conclusions obtained for the ID DM naturally hold 
for the 2D DM. We recall that the main difference between 
the ID and 2D geometry is that i) the spatial distribution 
of sources, which was constant in ID, is now q{r); and ii) 
the Galaxy has a side-boundary at a radius taken to be 
i? = 20 kpc. As a check, wc first used the 2D solution 
(presented in Appendix IA.2P with R = 20 kpc, but set q{r) 
to be constant. The best-fit parameters were in agreement 
with those obtained from the ID solution. We present in 
Table [5] the best-fit parameters for models II and III for 
L = 4 kpc in the 2D s olution where qjr) fol l ows th e SN 
remnant distribution of ICase fc Bhattacharval ([l99l). The 
values for the ID solution are also reported for the sake 
of comparison. The main difference is in the value of Kq, 
which varies by ~ 10% and also affects Va (by means of the 
ratio Va/\/Ko, which is left unaffected). This is c onsistent 
with the variations found bv lMaurin et al.l (I2002D . 



5. Results for radioactive species (free halo size L) 

We now attempt to lift the degeneracy between the halo 
size and the normalisation of the diffusion coefficient, using 
radioactive nuclei. The questions that we wish to address 
are the following: i) With existing data, how large are the 
uncertainties in L for a given model? ii) Do radioactive nu- 
clei provide different answers for models with different S? 
iii) Is the mean value (and uncertainty) for L obtained from 
a given isotopic/elemental ratio consistent with or stronger 
constrained than that obtained from another measured iso- 
topic/elemental ratio? iv) How does the presence of a local 
underdense bubble (modelled as a hole of radius rh, see 
Sect. 12. 3p affect the conclusions? 

Until now, almost all studies have focused on the 
isotopic ratios of ^^Be/^Be, ^^Al/^' ^Al ^eg/q^ ^nd 
^^Mn/Mn. An alternative, discussed in lWebber fc Soutoull 
(|1998I) . is to consider the Be/B, Al/Mg, Cl/Ar, and Mn/Fe 
ratios. The advantage of considering these elemental ratios 
is that they are easier to measure than isotopic ratios, and 
thus provide a wider energy range to which we can fit the 
data. Taking ratios such as Be/B maximises the effect of 
radioactive decay, since the numerator represents the de- 
caying nucleus and the denominator the decayed nucleus. 
However, the radioactive contribution is only a fraction of 
the elemental flux, and HEAO-3 data were found to be less 



Table 7. Most probable values for models II and III for 

the free parameters of the local bubble radius and/or 
the Galactic halo size L (constrained by B/C and ^°Be/^Be 
data). 



A'o X 10" S Va Va I ~ 

(kpc"Myr~^) (kms~^) (kms~^) (kpc) (pc) 



II 


8.6t°-J 


Q 2oq+o.oo5 






[4] 




II 




Q 2S4+0 006 






5 2+°-"^ 

>^-^-0.6 




II 




n nqr+O.OOS 

U.ZOO_g QQ4 








q+70 


III 




0.86l«:« 




351^ 


[4] 




III 


6 


0.86tH? 


19.4l«:^ 


13511J 


46t^ 




III 


o.stS., 


0.86iH! 




55^21 


8±? 


120±^EJ 



constr aining that the isotopic ratios in lWebber &: Soutoull 
(IT998I) . 

Below, we consider and compare the constraints from 
both the isotopic ratios and the elemental ratios. The 
data used are described in Appendix ID. 21 We discard 
^^Mn because it suffers more uncertainties than the oth- 
ers in the calculation (and also experimentally) due to the 
electron capture decay channel. The free parameters for 
which we seek the PDF are the four transport parameters 
{i^o, (5, Vc, Va\, plus one {L\ or two geometrical param- 
eters {L, r/j}, depending on the configuration considered. 
The main results of this section are thus in identifying the 
PDF of L for the standard DM, and the PDFs of both L 
and r/i for the modified DM. 

5.1. PDFs of L and using isotopic measurements 

Wc start with a simultaneous fit to B/C and ^°Be/^Be, for 
both model III (diffusion/convection/reacceleration), and 
model II (diffusion/reacceleration), the latter being fre- 
quently used in the literature. 

5.1.1. Simultaneous fit to B/C and ^°Be/^Be 

The marginalised posterior PDFs of L and rh and the cor- 
relations between these new free parameters and the prop- 
agation parameters of models II and III are given in the 
Figs, inland [71 respectively. The most probable values of the 
parameters are gathered in Table jT] 

For all configurations, the diffusion slope S and the 
Galactic wind Vc are unaffected by the addition of the free 
parameters L and rh- The B/C fit is degenerate in Kq/L 
and Va/VK, so that the values of Kq and Va vary as L 
varies. For model III, their evolution follows the relations 
given in Fig. [5] This implies that there is a positive cor- 
relation between Kq and Va, and Kq and L, as seen from 
Figs, m and |7l The uncertainty in the diffusive halo size L 
is smaller for II than for III. This is a consequence of the 
inclusion of the constant wind, which decreases the resolu- 
tion on Kq from 2% (Model II) to 10% (Model HI)— see 
e.g.. Table [5] or [7] — hence broadening the distribution of L. 

Below, the results for the standard DM — for which 
is set to be — and those for the modified DM — for which 
rfi is left as an additional free parameter — are discussed 
separately. This allows us to emphasise the impact of 



9 



A. Putze et al.: An MCMC technique to sample transport and source parameters of Galactic cosmic rays. II. 




Model II 




0.05 0.1 0.15 



0.05 0.1 0.15 



Fig. 6. Model II (diffusion/reacceleratioii): marginalised posterior PDF of the diffusive halo size L (right panels of the first 
and second row) and the local bubble radius (right panel of the last row) for the standard (r/i = 0, first row) and the 
modified DM (r/j ^ 0, second and last row), as constrained by using B/C and ^"Be/^Be data. The correlations between 
the geometrical parameters L and and the transport parameters 5, Kq, and Va, are shown in the 2D histograms. The 
colour code corresponds to the regions of increasing probability (from paler to darker shades), and the two contours 
(smoothed) delimit regions containing respectively 68% and 95% (inner and outer contour) of the PDF. 
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on the other parameters, which is different for models II Standard DM {r^ — 0): the parameter L is constrained to 
and III. be between 4.6 and 5.9 kpc for model II, having a most 

probable value at 5.2 kpc — a result compatible with other 
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studies (jMoskalenko et all [200lh — the posterior PDF of L 
extends from 25 to 85kpc for model III (most probable 
value at 46kpc). In terms of statistics, the best-fit model is 
stil model III, for which the x^/d.o.f is 1.41. 



Modified DM [rh 7^ 0): the presence of a local bubble re- 
sults in an exponential attenuation of the local radioactive 
flux, see Sect. 12.31 and Eq. ([TU)) . We thus expect to have 
a different best-fit parameter for L in that case. The re- 
sulting posterior PDFs of L and rh and the correlations 
to the propagation parameters for this modified DM are 
given in Figs. IH] and [7] for models II and III respectively. 
The most probable values are gathered in Table [7] (third 
and last lines). 

As expected, the local bubble radius rh is negatively 
correlated with the Galactic halo size L. The effect is more 
striking for model III, where the favoured range for L ex- 
tends from 1 to 50kpc. The most probable value is Sj^^kpc 

for a local bubble radius rh = 120j^2oPC- The x^/d.o.f of 
this configuration is 1.28, instead of 1.41 for the standard 
DM. The improvement to the fit is statistically significant 
according to the Fisher criterion. 

The situation for model II is different. The halo size 
L is already small for the standard configuration rh = 0. 
Adding the local bubble radius to the fit decreases the 
most probable value of L only slightly to 4^^ kpc and the 
measured value of rh is compatible with pc. In addition, 
the x^/d.o.f is 3.69 and hence poorer than for the configu- 
ration without the local bubble feature. In this model (dif- 
fusion/reacceleration, no convection), a local underdensity 
is not supported. 



5.1.2. Results and comparison with fits to ^^AI/^''AI and 
36CI/CI 

We repeat the analysis for the remaining isotopic ratios. 
The resulting marginalised posterior PDFs of the Galactic 
halo size L and the local underdensity rh are given in Figs.jS] 
and [9] for models II and III, respectively. The correlation 
plots with the transport parameters are similar to those of 
Figs. Eland [7] and are not repeated. 



Standard DM {rn = 0): as for the ^°Be/^Be ratio (red- 
dotted line) , L is well constrained in model II at small val- 
ues for the ^^A1/^''A1 (green-long dashed-dotted line), and 
^^Cl/Cl (blue dashed-dotted line) ratios, covering slightly 
different but consistent ranges from 4 to 14 kpc. The width 
of the estimated PDFs increases when moving from the 
^°Be/^Be ratio to the '^^Cl/Cl ratio, due to the decreasing 
accuracy of the data. In the same way, the adjustment to 
the data becomes poorer, as expressed by the increase in 
X^/d.o.f from 3.59 to 4.09. Used alone, the radioactive ratio 
^°Be/^Be constrains the most precisely the halo size L, but 
the constraints obtained with the other radioactive ratios 
arc completely compatible within the 2cr range. The most 
likely value of L is ascertained when the three radioactive 
ratios are fitted simultaneously (black solid line). 

The best-fit model is model III, where the overall cov- 
ered halo size range extends from 20 to 140 kpc. The most 
probable value found for L with 68% confidence level (CL) 
errors is 62]']}q kpc. 
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Fig. 8. Model II: marginalised posterior PDFs of the 
Galactic geometry parameters for the standard DM [rh — 
0, top panel) and for the modified DM [rh 7^ 0, bot- 
tom panels). The four curves result from the combined 
analysis of B/C plus isotopic ratios of radioactive species: 
B/C-^iOBe/^Bc (red dotted line), B/C+^SAl/^^Al (green 
long dashed-dotted line), and B/C-|-^^C1/C1 (blue dashed- 
dotted line) . The black solid curve represents the extracted 
PDF resulting from a simultaneous fit of B/C plus all three 
isotopic ratios. All PDFs are smoothed. 
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Fig. 9. Same as in Fig. |51 but for model III. 



Modified DM {rh ^ 0): the resulting marginalised poste- 
rior PDFs of L are shown in Figs. [5] and [3] (lower left) for 
models II and III, respectively. Again, the extracted PDFs 
for all radioactive ratios are completely compatible for both 
models. As described above, the decrease in L is more pro- 
nounced for model III than for model II. This decrease can 
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be observed for all radioactive ratios, independently of the 
model chosen. 

The resulting marginalised posterior PDFs of rh are 
given in Figs. [5] and IHl (lower right) for models II and III, 
respectively. The addition of an underdensity in the local 
interstellar medium is preferred by the data in the best-fit 
model III, but it is disfavoured in model II. The most prob- 
able values for Vh range from 90 pc for the '^^Cl/Cl ratio to 
140 pc for the the ^^A1/^''A1 ratio, and the overall fit points 
to a most probable radius of 130^20 P*^- 

These resul ts confirm an d exte nd the slightly differ- 
ent analysis of iDonato et al.l ()2002|) , who found that for 
model III, the best-fit values for was ^ 80 pc (see also 
Appendix [B]) . 

5.1.3. Envelopes of 68% CL 

Confidence contours (for any combination of the CR fluxes) 
corresponding to given confidence levels (CL) in the dis- 
tribution can be drawn, as detailed in Appendix A and 
Sect. 5.1.4 of Paper I. From the MCMC calculation based 
on the B/C + lOBe/^Be + 26^1/27^1 + 36^/^1 constraint, 
we select all sets of parameters for which the meets the 
68% confidence level criterion. For each set of these pa- 
rameters, we calculate the B/C and the three isotopic ra- 
tios. We store for each energy the minimum and maximum 
value of the ratio. The corresponding contours (along with 
the best-fit ratio) for models II (standard DM, red) and III 
(standard and modified DM, blue) are drawn in Fig.[TUl To 
ease the comparison with the data, all results correspond 
to IS quantities (the approximation made in the demodu- 
lation procedure, see Paper I, is negligible with respect to 
the experimental error bars). 

We see that the present data already constrain very well 
the various ratios for the standard DM. The difference be- 
tween the results of models II and III arc more pronoimccd 
at high energy (effect of (5), as seen from the B/C ratio be- 
yond 10 GeV/n. All contours are pinched around 10 GeV/n, 
which is a consequence of the energy chosen to renormalisc 
the flux to the data in the propagation code. In princi- 
ple, the source abundance of each species may be set as 
an additional free parameter in the fit (Paper I), but at 
the cost of the computing time. The three isotopic ratios 
(lOBe/^Be^e, Al/^^Al, and ^^Cl/Cl) provide a fair match to 
the data for all models, considering the large scatter and 
possible inconsistencies between the results quoted by vari- 
ous experiments. In particular, for ^°Be/^Be, new data are 
necessary to confirm the high value of the ratio measured 
at ~ GeV/n energy. 

The envelope for the modified DM is quite large at 
high energy, because the uncertainty in is responsi- 
ble for a larger scatter in the other parameters. The two 
standard DM contain non-overlapping envelopes beyond 
GeV/n energies. This means that to disentangle the mod- 
els, having measurements of the above isotopic ratios in 
the 1 — 10 GeV/n may be more important than just having 
more and higher quality data at low energy. 

5.1.4. General dependence of L with 6 (for rh = 0) 

To investigate the difference in the results obtained from 
models II and III, we fit B/C and the three isotopic ratios 
for different values of 5 (a similar trend with L is obtained 



if just one isotopic ratio is selected). The analysis relies on 
the Minuit minimisation routine to quickly find the best-fit 
values, as described in Maurin et al. (2 01Q) . The evolution 
of the parameters with 5 is shown on the left side of Fig.fTTl 
The bottom panel shows the evolution of xim^/f^-O-i, where 
we recover that the best-fit 5 for model II (dashed-blue 
line) lies around 5 k, 0.2, whereas that for model III (solid- 
black line) lies around 5 ~ 0.8. As already underlined, the 
contribution to the value is dominated by the B/C con- 
tribution because as discussed in Appendix [Cj the values 
of transport parameters that reproduce the B/C ratio are 
expected to remain within a narrow range. This explains 
what is observed in the various panels showing these com- 
binations. For model II, we emphasise that for 5 > 0.5, 
the best-fit value for Va is zero (Model II becomes a pure 
diffusion model). 

The most important result, given in the top panel, is for 
L as a function of 6, where any uncertainty in the determi- 
nation of 5 translates into an uncertainty in the determina- 
tion of L. When a Galactic wind is considered (Model III, 
black-solid line), the correlation between L and S is stronger 
than for model II (no wind). There is no straightforward 
explanation of this dependence. The flux of the radioactive 
isotope can be sho wn to be A^''^'^(0) sa hq/^/Kyr^ (e.g., 
|m aurin et all l2006l ) . Since secondary fluxes should match 
the data regardless of the value for 5, this implies that the 
ratio ^°Be/^Be depends only on y/ Kjtq . At the same time, 
to ensure that the secondary-to-primary ratio is constant, 
we must maintain a constant L/K. The difficulty is that the 
former quantity is a constant at low rigidity where the iso- 
topic ratio is measured, whereas the latter quantity should 
remain as close a possible to the B/C data over the whole 
energy range. Hence, all we can say is that the variation 
in L with S is related to the variation in Kq/L with S, as 
shown in the second figure (left panel) of Fig. [TT] 

We note th at all the calcu l ations in the paper are based 
on the W03 (jWebber et al.l l2003l) fragmentation cross- 
sections. The impact of using the W03 set or the GAL09 
set (provided in the widely used GALPROP packagcQ) on 
the determination of the halo size L is shown as thin-dotted 
lines (left panel, same figure)0- Any difference existing be- 
tween these two sets of production cross-sections has no 
impact on the best-fit value for L: thin-dashed curves (ob- 
tained with GAL09 cross-sections) almost match the thick- 
solid curves (obtained with W03 cross-sections). For other 
ratios, the effect of the GAL09 cross-sections is always the 
same, so it is not discussed further. 

5.1.5. General dependence of L with 5 (for r/, 7^ 0) 

We repeat the analysis with the underdensity as an ad- 
ditional free parameter. The dependence of L and on 
the diffusion slope 5 is shown in the right panel of Fig. [TT1 
A comparison between the left and the right panel shows 
that the combinations of parameters Kq/L, Vaj \JKq^ and 
Vc are almost unaffected by the presence of a local bubble; 
the x^/d.o.f is also only slightly affected. 

For (5 < 0.2, r/j is consistent with for both 
model II (diffusion/reacceleration) and model III (diffu- 

^ http://galprop.stanford.edu/wcb_galprop/galprop_homc.htnil 
^ The impact o n the tr ansport parameters is detailed in Sect. 7 

of iMaurin et all (1201^ ): the region of the best-fit values is 

slightly displaced, as seen in the figure. 
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Fig. 10. Shown are the envelopes of 68% CL (shaded areas) and best-fit (thick lines) ratios for the standard DM II 
{rh = 0, red) and for model III (standard and modified DM, blue) in the ID geometry (based on the B/C + ^°Be/^Be 
^^A1/^''A1 + ■^^Cl/Cl constraint). All quantities arc IS. The data arc demodulated using the approximate procedure 
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sion/convcction/reacceleration). For model II, the size of 
r/i suddenly jumps to ^ 100 pc. But for S > 0.3, it returns 
to the pure diffusion regime, decreasing abruptly (to 
a non- vanishing value) and L becoming vanishingly small. 
In this regime, the thin-disc approximation is no longer 
valid and nothing can be said about it. For model III, the 
plateau ^ 100 pc is stable for all 6 > 0.2. The under- 
dense bubble also stabilises the value of the halo size L. 
The way of understanding this trend is as for the stan- 
dard DM, but now the flux of the radioactive species reads 
A^^',^'^(0) oc exp(-r/i/VA'7To)/Vif7ro. The weaker depen- 
dence of L with S must be represented by this formula. 
We underline that for all best-fit configurations leading to 
rii 0, the improvement is statistically meaningful com- 
pared to the case r/j = 0. 



5.2. Isotopic versus elemental measurements 



A similar analysis can be carried out using elemental ra- 
tios instead of isotopic ones. As before, the best-fit values 
of well-chosen combinations of the transport parameters 
{Kq, S, Vc, Va} are left unchanged when radioactive species 
are added to the fit (same values as in Fig. [TT|) . 
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Fig. 12. Best-fit value of the halo size L as a function of 5 
in standard DM, based on a fit on B/C plus a ratio where 
a radioactive species is present: B/C-f-Bc/B (black small 
symbols), B/C-t- Al/Mg (blue medium-size symbols), and 
B/C-f Cl/Ar (pink large symbols). The dashed lines (square 
symbols) refer to model II, and the solid lines (circles) refer 
to model III. 



5.2.1. General dependence of L with 6 

For the standard DM {rh = 0), the dependence of the diffu- 
sive halo size L on the diffusion slope S is shown in Fig. [T^ 
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Fig. 11. Left panel: standard DM model {rh = 0) — thin-dotted lines are derived using the GAL09 instead of the W03 
fragmentation cross-sections. Right panel: modified DM model {r^ ^ 0). For both panels, shown are the best-fit parame- 
ters on B/C -I- i°Be/^Be -I- 2''A1/2''A1 -I- ^'^Cl/Cl data, as a function of the difl:usion slope b. The latter is varied between 
0.1 and 1.0 for model II (blue lines, open and filled squares) and model III (black lines, open and filled circles). From 
top to bottom, L. Kq/L. Va/\/KQ, and Vc as a function of 6 are shown. The bottom panel shows the best x^/d.o.f for 
each 5. 



for the three combinations B/C + Bc/B, B/C + Al/Mg, 
and B/C -I- Cl/Ar. The trend is similar to that for isotopic 
ratios: L increases with increasing 5. The main difference is 
that the increase is sharper for both models II and III. For 
the former, only a small region around 5 « 0.2 corresponds 
to small halo sizes. For the latter, the halo size increases 
sharply above 5 > 0.6. 



For completeness, similar fits were carried out for the 
modified DM {rh ^ 0). However, adding an additional de- 
gree of freedom only worsens the situation, and the models 
converge to arbitrarily small or high values of L and r/j. 
Finally, if we fit the combined B/C data, the three isotopic 
ratios and the three elemental ratios, we do not obtain more 
constraints than when fitting B/C and the three isotopic ra- 
tios. This may indicate that the models have difficulties in 
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Fig. 13. Same as in Fig.fTUlbut for the ratios Be/B, Al/Mg, 
and Cl/Ar. 



fitting all these data together: either the model is incom- 
plete or the data themselves may show some inconsisten- 
cies. This is more clearly seen from the comparison of the 
model calculation and the data for these elemental ratios 
(see below). 



the Be/B ratio and the two other ratios is that Be and B 
are pure secondary species, whereas all other elements may 
contain some primary contribution which can be adjusted 
to more closely match the data. This also explains why the 
Be/B ratio reaches an asymptotical value at high energy 
(related to the respective production cross-sections of Be 
and B), whereas the two others exhibit more complicated 
patterns. The low-energy Be/B ratio is related to cither the 
model or the energy biases in the production cross-sections 
for th ese elements (which is still possible, e.g. lWebber et all 
|2003[ ). or to systematics in the data. To solve this issue, 
better data over the whole energy range are required. 

5.3. Summary and generalisation to the 2D geometry 

Using radioactive nuclei in the ID geometry, we found 
that in model II (difFusion/reacceleration), L ^ 4 kpc 
and rfi ^ 0, and for the best-fit model III (diffu- 
sion/convection/reacceleration), L ^ 8 kpc and ~ 
120 pc. The halo size is an increasing function of the dif- 
fusion slope 6, but in model III the best-fit value for rh 
remains ~ 100 pc for any 5 > 0.3. This value agrees 
with direct observation of the LISM (see Appendix [B|) . 
Measurement of elemental ratios of radioactive species are 
not yet precise enough to provide valuable constraints. 

For now, there are too large uncertainties and too many 
inconsistencies between the data themselves to enable us to 
point unambiguously toward a given model. Moreover, one 
has to keep in mind that any best-fit model is relative to a 
given set of data chosen for the fit (see Sect. 14. 2|) . We note 
that there may be ways out of reconciling the low-energy 
calculation of the Be/B ratio with present data, e.g., by 
changing the low-en ergy form of the diffusion coefficient 

aurm et al.ll2010[) . but this goes beyond the goal of this 
paper. 

All these trends are found for the models with 2D ge- 
ometry. We calculate in Table [5] the best- fit parameters for 
the standard model II (r/i = 0) and the modified model III 
{^h 7^ 0). The values for the ID geometry are also reported 
for the sake of comparison. Apart from a few tens of percent 
difference in some parameters, as emphasised in Sect. 14.51 
some differences are expected if the size of the diffusive halo 
L is larger than the distance to the side boundary R, which 
is du = 12 kpc in the 2D geometry. It is a well-known 
result that the closest boundary limits t he effective diffu- 
sion r egion from where CR can originate (|Taillet fc MaurinI 
|2003| ). For model II, L is smaller than du. We obtain a 
smaller than 10% difference for Kq, and a ^ 30% difference 
for L. For the modified model III, the halo size has a larger 
scatter (see previous sections), with L^'lf* = 13.6 > d^. The 
geometry is thus expected to affect the determination of L. 
We find that ijif* = 4.3 and that the value of Kq is thus 

i5'i5*/-^2lf* ~ 3 times larger, and Va is ~ a/S times larger 
than in ID. 



5.2.2. Envelopes of 68% CL 

From the same set of constraint as in Sect. 15.131 fi.e.. B/C 
and the isotopic ratios of radioactive species only) , wc draw 
the CL for the elemental ratios in Fig. [TS] 

Given their large error bars, the elemental ratios are in 
overall agreement with the data, except at low energy and 
especially for the Be/B ratio. The main difference between 



6. Conclusions 

We have used a Markov Chain Monte Carlo technique to 
extract the posterior distribution functions of the free pa- 
rameters of a propagation model. Taking advantage of its 
sound statistical properties, we have derived the confidence 
intervals (as well as confidence contours) of the models for 
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Table 8. Best- fit parameters on 

B/C+i"Bc/SBe26+Al/2'^Al +^^C\/Ch ID vs 2D DM. 



Config. Ko X 10- 


S 




Va 




xVd.o.f 


(kpc^Myr^^ 


) 


(km 




(kpc) 




ID-II-L 15.4 


0.23 




92.5 


6.2 ■•■ 


3.09 


2D-II-L 14.9 


0.24 




90.8 


8.8 ■•■ 


3.04 


ID-III-Lrh 1.90 


0.83 


18.9 


73.5 


13.6 0.13 


1.43 


2D-III-Lrh 5.24 


0.85 


18.3 


123. 


4.3 0.16 


1.48 



fltrxes and other quantities derived from tire propagation 
parameters. 

In the first paper of this series (Paper I), we focused on 
the phenomenologically well-understood LBM to ease the 
implementation of the MCMC. In contrast, here we have 
analysed a more realistic DM. In agreement with previous 
studies, when B/C only is considered, we have confirmed 
that a model with diffusion/convection/reacceleration is 
more likely than the diffusion/reacceleration case. The for- 
mer would imply that 5 ^ 0.8, whereas the latter would 
imply that 6 ^ 0.2. This result does irot depend on the 
halo size: we provided simple parameterisations to obtain 
the value of the transport parameters for any halo size L. 
If mere eye inspection of the published AMS-01 data shows 
consistency with the HEAO-3 data (covering the same en- 
ergy region), a B/C analysis based on AMS-01 data (instead 
of IIEAO-3) also indicates that convection and reaccelera- 
tion is required, but now providing a diffusion slope d ^ 0.5, 
closer to theoretical expectations. Data from PAMELA or 
high-energy data from CREAM and TRACER are required 
to help solving the long-standing uncertainty in the value 
for S. 

A second important topic of this paper has been the 
halo size L of the Galaxy and the impact of the underdense 
medium in the solar neighbourhood. The determination of 
L is for instance crucial to predictions of antimatter fluxes 
from dark matter annihilations. The size of the local un- 
derdense medium is as important, as it can bias the de- 
termination of L. We provided a step- by-step study of the 
various radioactive clocks at our disposal. Our detailed ap- 
proach can serve as a guideline as how to take advantage of 
future high-precision measurements that will soon become 
available (e.g., from AMS). The main conclusions about the 
constraints provided by the radioactive species are, in diffu- 
sion/reacceleration models, L ~ 4 kpc and no underdense 
local bubble is necessary to match the data. For the best-fit 
model, which requires diffusion/convection/reacceleration, 
L ^ 8 kpc with rh ~ 120 pc. For both models, the halo 
size found is an increasing function of the diffusion slope 
S. A striking feature is that in models with convection, the 
best-fit value for remains ^ 100 pc for any S > 0.3. 
For instance, the B/C AMS-01 data (which implies that 
6 ~ 0.5) and the radioactive ratios are consistent with a 
wind and a local underdense bubble. This very value of 
rh ~ 100 pc is also supported by direct observation of the 
LISM (see Appendix 

As emphasised in this study, the determination of the 
value of L and r/j strongly depends on the value of S. 
For all these parameters, high-energy data of secondary-to- 
primary ratios, data in the ^ 1 GeV/n— 10 GeV/n range 
for isotopic ratios (of radioactive species), and/or data for 



the radioactive elemental ratios in the 1 — 100 GeV/n en- 
ergy range are necessary. This is within reach of several hy- 
ing and forthcoming balloon-borne projects and satellites 
(PAMELA, AMS). 
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Appendix A: Solutions of the diffusion equation 

We provide below the solutions for the diffusion equation 
with a constant wind Vc and a single diffusion coefficient 
K(E ) in the whole Ga laxy. In the ID version of the model 
(e.g.. lJones et al.ll200l|) . the source distribution and the gas 
density do not depend on r, so that the propagated fiuxes 
depend only on z. 

The derivation of these solutions is very similar 
an d has no add i tional difficulties to those experienced 
by iMaurin et al.l ()2001|) , to which we refer the reader for 
more details. As both frameworks (ID and 2D) exhibit 
similar forms, formulae are written for the ID model only. 
Formulae for the 2D case are obtained by replacing some 
ID quantities by their 2D counterparts, as specified below. 

A.l. ID-model 

The starting point is the transport equation 
(|Berezinskh et al.l Il990t) . We assume that the diffusion 
coefficient K does not depend on spatial coordinates. A 
constant wind Vc blows the particle away from the Galactic 
disc, along the z directio n. In the thin-disc approximation 
(e.g.. IWebber et anil992f ). the diffusion/convection for the 
ID-model (discarding energy redistributions) is 

{"^S + Vc-^ + r,ad + 2/irtot<5(z)| N = CN^ Qiz). 

In this equation, N is the differential density of a given 
CR species, Fj-ad = 1/(7''') is its decay rate, and Ftot = 
X^isM "iSMi'CisM is its destruction rate in the thin gaseous 
disc (nisM =H, He). The right-hand side (r.h.s) of the equa- 
tion is a generic source term, that contains one of the fol- 
lowing three contributions, i.e., Q{z) — 'P{z)+S{z)+TZ{z): 

i) Viz) = 2h6{z) X QqQIE) is the standard primary source 

term for sources located in the thin disc. The quantity 
Qq is the source abundance of nucleus j whose source 
spectrum is Q{E) cx ^''7^-". 

ii) S{z) = 2hS{z) X TP;^'NP{z = 0,E) is the standard sec- 
ondary source term (also in the disc), where T^^" = 
nvfjP^^ is the production cross-section of nucleus p into 
s. This simple form originates from the straight-ahead 
approximation used when dealing with nuclei (see, e.g., 
IMaurin et al.ll2001l for more details). 

in) TZ{z) = V^^^l^N^{z, E) described a contribution from a 
radioactive nucleus r, decaying into s' in both the disc 
and the halo. 

The equation is even in z so that it is enough to solve it 
in the upper-half plane. The use of the standard boundary 
condition N{z = L) = and continuity of the density and 
the current at the disc crossing completely characterises the 
solution. 
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A. 1.1. Stable species 



and 



For a mixed species, primary and secondary standard 
sources add up, so that, for a nucleus k with no radioactive 
contribution, the source term is rewritten as 

QZciE) = q'^Q{E) + rfot"^^(^ = 0, E), (A.2) 

and the corresponding equation to solve is then 

C^N^ = 2h5{z) ■ Q2.,{E) . 

We find the solution in the halo, apply the boundary con- 
dition N[z = L) = 0, and then ensure continuity between 
the disc and the halo, so that 

7V™(z) = 7V™(0) • exp(^=^/2^) si"M^"'(^- ^)/2) ^^^^ 
^ ' y I sinh(S'™L/2) ^ ' 



and 



7V"(0) = 



2feQdisc(^) 



The quantities S™ and A™ are defined as 



^" = K + +/\5"coth 



(A.4) 

(A.5) 
(A.6) 



disc 



cosh(5JL/2) 



(A.9) 



- sinh(S"'i/2) 
cosh(5''L/2) 



a K(2 - 



2/7P 



a^ ( ^c(2-— )+2/ir-' 



where 



,kVcL/2 



A; 



A 



2i5:'' 



rad 



rad 



K3 



The superscript on A' indicates that the diffusion coefficient 
is to be evaluated at a rigidity calculated for the nucleus m. 
The latter can differ from one nucleus to another because, 
the calculation is performed at the same kinetic energy per 
nucleon for all the nuclei (hence at slightly different rigidi- 
ties for different nuclei) . To compare with the data, the flux 
is calculated at z = 



N\{)) = m + Q 



X sinh( ^— ^ ) -I- A cosh( ^— ^ ) - 



cosh ^ 



(A.IO) 



A. 1.2. Adding a /3-decay source term: general solution 

It is emphasised in lMaurin et al.l (|200l[ ) that the i''Be-J>i°B 
channel contributes up to 10% in the secondary boron flux 
at low energy and cannot be neglected. Although the spa- 
tial distribution of a radioactive nucleus decreases expo- 
nentially with z, we have to consider that the source term 
is emitted from the halo, complicating the solution. The 
equation to solve for the nucleus j, which is /3-fed by its 
radioactive parent r is 



where N^{z, E) is given by Eq. (|A.3p . The solution is found 
following the same steps as above, although it has a more 
complicated form (due to a non-vanishing source term in 
the halo). 

If we take into account both the standard source term 
*9disc(-^) ^'^'i radioactive contribution of the nucleus 
A^'', we obtain: 



N^{z) 



sinh(S'^'(i-z)/2) 
smh{S3L/2) 
_ cosh(g^z/2) I 

cosh(5JL/2) J ^ 

e <! Asinh ( +Acosh(y(L-z) 



X exp 



(V,z/2K-') 



where 



e 



rad 



Ar'-(O) 



A:j(A2 - A2) sinh(S''-L/2) 



(A.7) 



(A.8) 



A. 1.3. Solution including energy redistributions 

When energy redistributions are included, the solution 
N^{z) in the halo remains the same because our model 
assumes no energy redistributions in that region. Only the 
last step of the calculation changes (ensuring continuity 
during the disc crossing). The new solution is denoted A/'(0). 

For the case of a mixed species m without radioactive 
contribution, the result is straightforward: the solution for 
the halo is still given by Eq. (jA.3|) . but A/'™(0) is now given 
by 



AA™(0)=Af™(0)- 



2h 



dE 



dE^ 



which is solved numerically, ^^'"(0) being the solution when 
energy terms are discarded, i.e., Eq. (jA.4|) . The terms a{E) 
and b{E) describing energy losses and gains are discussed 
in Sect. HHl 

When a radioactive contribution exists, the constant left 
to determine is zu from Eq. (jA.7|) . which we denote now 

-* = --^(6(^^)^ + c(£;)^^j . (A.11) 

As above, w denotes the quantity evaluated without energy 
redistribution, whereas Af-' (0) denotes the equilibrium flux 
at 2 = 0. To ensure J\f^{0) also appears in the l.h.s. of the 
equation, we form the quantity 



s = e 



A sinh( ^^-^ ) + A cosh( ^^-^ ) - 



i/A 



cosh ^ 



(A.12) 
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Hence N"^{0) = m + E, and we can add to both sides of 
Eq. (jA.lip the quantity S, so that we recover the standard 
form 

^'<»)^-<°)-§('<-)^-<-)^)^ 

which we solve numerically. 

This is the solution in the disc {z = 0). The solution for 
any z is obtained from Eq. (|A.7p . making the substitution 

■uu ^ U7* ^ Af^O) ~ E. (A.13) 

We note that for 8 = (i.e., no radioactive contribution) 
the result for standard sources in the disc is recovered. 

A. 2. 2D geometry 

Cylindrical symmetry is now assumed, both the CR den- 
sity N and the source terms depending on r. Compared to 
Eq. (jA.ip . the operator now acts on N{r, z). 

An expansion along the first order Besscl function is 
performed 

oo 

Ar(r,z) = ^7V,(z)Jo(o-^) . (A.14) 

i=l 

The quantity Q is the i-th zero of Jq, and this form auto- 
matically ensures the boundary condition N{r = R,z) ^ 0. 
We have 

-A.Jo (C4) = ^ Jo (C4) , 

so that each Bessel coefficient Ni{z) follows an equation 
very similar to Eq. (jA.l|) . where 

Trad Ejad 4* ' 

and where each source term must also be ex panded on the 
Bessel basis. More details can be found in iMaurin et al.l 

(pool . 

The full solutions for mixed species, with stable or ra- 
dioactive parents, is straightforwardly obtained from ID 
ones, after making the substitutions 

7VJ"(z)™^''^'7V/(z), (A.15) 

^l^\/^+4f +4^ , (A.16) 

A= 2D^dci ^ ^ + K + KSj coth 1^^^ , (A.17) 
and 

e^{s\N^m^ei{s:,N:{o)) , (A.18) 

u,^S^,A^)=^uji{SUAi), (A.19) 

A'-(5'')=^A[(5[) . (A.20) 

The above formulae, for the ra dioactive source, dif - 
fer slightly from those presented in iMaurin et all (1200111 . 
However, the only difference is in the flux for z 7^ 0, which 
was not considered in this paper. 



Appendix B: The local bubble 

The underdensity in the local interstellar matter (LISM) 
is coined the local bubbl43- The LISM is a region of ex- 
tremely hot gas 10^ — 10^ K) and low density (n < 
0.005 cm~'^) within an asymmetric bubble of radius < 
65-250 pc s urrou n ded by dens e neut r al hydrogen walls 
(iSfeir et al.l [1991 iLinskv et al.l [2000I: iRedfield fc Lhlslol 
|200G[ ). This pict ure has been r efined by subsequent stud- 
ies, e.g.. LaUeme nt et al.l (|2003f ). The Sun is located inside 
a local interstellar cloud (LIC) of t ypical extension ~ 50 
pc whose density Nhi f -^ 0.1 cm"^ (jCloeckler et al.ll200l 
iRedfield fc Falcoiil l2008l ) . Despite these successes, a com- 
plete mapping and understanding of the position and prop- 
erties of the gas/cloudlets filling the LISM, as well as the 
issue of interfaces with other bubbles remains cha llenging 
re.g.. lRedfield fc Linskvl[2008t iReis fc Corradill2008D . Based 
on existing data, numerical simulations of the local bubble 
infer that it is the result of 14 — 19 SNe occurring in a 
moving group, which passed th rough the present day local 
Hi ca vity 13.5 — 14.5 Myr ago (jBreitschwerdt fc de Aville^ 
I2006h . The same study suggests that the local bubble ex- 
panded into the Milky Way halo roughly 5 Myr ago. 

A last important point, is that of the existence of tur- 
bulence in the LISM to scatter off CRs. The impact of the 
underdense local bubble on the production of radioactive 
nuclei as modelled in Eq. (fTU)) depends whether the trans- 
port of the radioactive nuclei in this region is diffusive or 
not. In a study based on a measurement of the rad io scintil- 
lation of a pulsar located within the local bubble, ISpangleil 
(|2008l ) infers that values for the line of sight component of 
the magnetic field are only slightly less, or completely con- 
sistent with, lines of sight through the general interstellar 
medium; the turbulence is unexpectedly high in this region. 

These pieces of observational evidences support the 
model used in Sect. 12.31 leading to an enhanced decrease 
in the flux of radioactive species at low energy. A detailed 
study should take into account the exact morphology of the 
ISM (asymmetry, cloudlets). However, there are so many 
uncertainties in this distribution and the associated level of 
turbulence, that a crude description is enough to capture a 
possible effect in the CR data. 

Appendix C: MCMC optimisation 

The efficiency of the MCMC increases when the PDFs of 
the parameters are close to resembling Gaussians. Large 
tails in PDFs require more steps to be sampled correctly. 
A usual task in the MCMC machinery is to find some com- 
binations of parameters that ensure that these tails disap- 
pear. This was not discussed in the case of the LBM as the 
efficiency of the PDF calculation was satisfactory. In ID (or 
2D) DMs, the computing time is longer and the efficiency is 
found to be lower. To optimise and speed up the calculation, 
we provide combinations of parameters that correspond to 
a Gaussian distribution. 

A typical PDF determination with four free parame- 
ters {Vc, 5, Kq, Va} (see next section) is shown in Fig. IC.ll 
The diagonal of the left panel shows the PDF of these pa- 
rameters (black histogram), on which a Gaussian fit is su- 
perimposed (red line). We see a sizeable tail for the Kq 

* For a state-of-the-art view on the subject, the reader is re- 
ferred to the proceedings of a conference held in 2008: The Local 
Bubble and Beyond II — Jrttp://lbb. gsfc.nasa.gov/ 1 
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Binary space partitioning step 

1 Chain analysis ^ 



Burn-In length; 2.9 
Correlation Length: 5.8 
^ifd- 135156/400000 = 33.789% 
X', /d.o.f : 1.4 7 




Binary space partitioning step 

chain analysis ^ 



Burn-In length; 3.0 
Correlation Length: 1.9 
fi„j: 209947/400000 = 52.487* 
X' /d.o. f : 1.47 




0.8 0.9 



0.004 0.006 O.OC 



Fig. C.l. Posterior PDFs of the model parameters (using the Binary Space Partitioning step — see Paper I, and the B/C 
constraint). The diagonals show the ID marginalised PDF of the indicated parameters, and the red line results from a 
Gaussian fit to the histogram. Off-diagonal plots show the 2D marginalised posterior PDFs for the parameters in the 
same column and same line, respectively. The colour code corresponds to the regions of increasing probability (from paler 
to darker shades), and the two contours (smoothed) delimit regions containing respectively 68% and 95% (inner and 
outer contour) of the PDF. Left panel: PDFs for {Vc, S, Kq, Va}- Right panel: the same PDFs but shown for a different 
combination of the parameters {Vc, S, Kq/L x 50'^, K/ \/Ka x 3(5(4 — 6^){4: — 6)}. 



parameter, and a small asymmetry for the Va parameter. 
The right panel shows the same PDFs, but for the following 
combinations of the transport parameters: 



Kn 



Va 



1j 



X 3(5(4-52)(4_ J) 



(C.l) 



(C.2) 



These forms are inspired by the known degeneracies be- 
tween parameters. For instance, in diffusion models, the 
secondary to primary ratio is expected to remain un- 
changed as long as the effective grammage (a;) of the 
model is left unchanged. For pure d iffusion, we obtain 
(jJones et all l2001t iMaurin et all |2006[) for the grammage 

(x) = i:vL/{2K) = i:L/{2cKo(n/lGYf), where E is the 
surface density. Apart from the Kq — L degeneracy, the 
parameters Kq and S are correlated. We find that the com- 
bination Kq X 50* is appropriate for removing the Kq PDF's 
The origin of the value 50 is unclear. It 



tail (see Fig. IC.lj) 
may be related to the energy range covered by B/C HEAO- 
3 data on which the fits are based. The combination used 
for Va [Eq. (|C.2p ] comes directly from the form of the reac- 
celeration term Eq. ([S])- Rcacceleration only plays a role at 
low energy, so we can take TZ^ « 1 and end up with the 
combination presented in Eq. (jC.2p . The independent ac- 
ceptance /ind) defined in Paper I as the ratio of the number 
of independent samples to the total step number, increases 
from 1/3 to 1/2 by using the above described parameter 
combinations for the four parameter model presented in 

Fig, [cm 



A last combination is for the local bubble parameter : 

rii 



(C.3) 



This comes from the form of Eq. ([TU)) , where the flux damp- 
ing for radioactive species (due to the local bubble) is ef- 
fective only at low energy (7 « 1, TZ ~ 1). 

Appendix D: Datasets for CR measurements 

D.l. B/C ratio 

Unless specified otherwise, the reference B/C dataset used 
throughout the paper is denoted data set F: it consists of i) 
low-e nergy data taken by the IMP7-8 (iGarcia-Munoz et all 
I1987D . the Vovaee r lfc2 (iLukasiak et al.l I1999D . and 
the ACE-CRIS (|de Nolfo et al.l I2006D spacecrafts; 
i i) intermediate energ ies acquired by HEAO-3 data 
(lEngelmann et aT 19901] : and iii) higher energy data from 
Space lab (ISwordy et anil990[ ) and the pubhshed CREAM 



data ()Ahn et al. 2008 ). Other existing data are discarded 
either because of their too large error bars, or because of 
their inconsistency with the above data (see Paper I). 



D.2. Isotopic and elemental ratios of radioactive species 

For ^°B e/^Be, the d ata are taken from balloon 
flights (iHagen et al.l 119771: iBufflngton et al.l 119781: 
IWebber fc Kisbl 119791) including t he ISOMAX balloon- 



borne instrument 
the IMP-7/8 



(iHams et al.l 



.20041) ■ and from 

(IGarcia-Mun oz et al.l 119771) . ISEE-3 



(jWiedenbeck fc Greined 119801) . Ulvsses (IConneUl Il998l) . 
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Vovager (iLukasiak et all Il999t) . and ACE spacecrafts 
(jYanasak et al.ll200l[) . Foi ^^Al/^'^Al, th e data consist of 
a series of balloon flig hts (iWebbed [1981 . and the ISEE- 



3 (IWiedenbeckl Il98.tl. Vo yager_ (|Lukasiak et all I1994D . 

Conneli l Il998[) . and ACE spacecrafts 



Ulysses (ISimpson 



(jYanasak et al 



20011) . For 36C1/C 1, the data are from 

"aDQ ,^ 

and ACE (jYanasak et al.l 



the CRI SIS balloon (lYounget al.l ll98lD . and from the 
Ulyss es (jConneU et al.l I1998D 



l2001f ) spacecrafts 

The dat a for the elemen tal ratios come from the HEAO- 
3 (jE ngclm ann et al.l Fl990[ ). Ul ysses (jP uvcrn ois fc Thayer! 
[199^, and the ACE spacecrafts Ide Nolfo et alTl2006l). The 
publi shed ACE data on Al/Mg and Cl/Ar (jCeorge et aD 
l2009f ) were not used as Be/B is not provided. 
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